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SUMMARY 

The surface temperature distributions due to thermocaplllary convections 
In a thin liquid layer with heat fluxes Imposed on the free surface are Inves- 
tigated. The nondlmenslonal analysis predicts that, when convection Is 
Important, the characteristic length scale In the flow direction L, and the 
characteristic temperature difference AT 0 , can be represented by 
L ~ (A2Ma) 1 ^Lp and AT 0 ~ (A 2 Ma) _1 / 4 ATR, respectively, where Lp and ATp 
are the reference scales used In the conduction-dominant situations with A 
denoting the aspect ratio and Ma the Marangonl number. Having had L and 
AT 0 defined, the global surface-temperature gradient (AT 0 /L), the global 
thermocaplllary driving-force, and other Interesting features can then be 
readily determined. Finally, numerical calculations Involving a Gaussian heat 
flux distribution are presented to Justify these two relations. 


INTRODUCTION 

The availability of a low-gravity environment aboard a spacecraft opens 
new possibilities for the growth of better quality crystals. Among various 
techniques used for crystal growth, contalnerless processes such as the 
floating-zone method seem to have a very great potential to benefit from the 
space environment. 

In the floating-zone method, a nonuniform heat flux Imposed along the free 
surface of the melt may Induce a thermocaplllary flow In the melt If the sur- 
face Is reasonably free from contamination. Chang and Wilcox (1975, 1976) 
Investigated floating-zone problems numerically by assuming a vertical cylin- 
drical melt with planar solid-liquid Interfaces. Clark and Wilcox (1980) later 
pointed out and corrected certain errors In Chang and Wilcox's numerical formu- 
lation. The main feature of the numerical work Is the prediction of flow cells 
(toroidal vortices) that extend through the entire fluid. Unfortunately, these 
numerical studies were limited to cases with small Marangonl numbers (1 to 10) 
which are too small to simulate actual situations, which are typically on the 
order of magnitude of greater than 10 2 . 

In order to simplify the analysis, most of the studies In the past dealt 
with the flow phenomena In a half-zone configuration (Chun, 1980a - Lai, 1984; 
Ostrach, Kamotanl, and Lai, 1985; Prelsser, Scharmann, and Schwabe, 1981 - 
Schwabe and coworkers, 1978; Vargas, Ostrach, and Kamotanl, 1982). In these 
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models a liquid bridge Is suspended vertically between two circular rods with 
a higher temperature at the upper rod and a lower temperature at the bottom 
one. Several Investigators have also studied the flow phenomena In a two- 
dimensional slot with different temperatures specified on the side walls (Lai, 
1984; Sen and Davis, 1982; Stranl, Plva and Grazlanl, 1983). The thermocapll- 
lary convections In a two-dimensional thin liquid layer with specified surface 
temperature distributions were also studied In the past (Lai, 1984; Plmputka 
and Ostrach, 1980; Smith and Davis, 1983a, 1983b). However, the coupling among 
the flow field, the temperature field and the Imposed heat flux was not taken 
Into consideration. In order to demonstrate certain features due to such a 
coupling, a simple model Is being proposed In this Investigation. 

The proposed model Is essentially a two-dimensional thin liquid layer with 
a nonuniform heat flux Imposed on the free surface. It differs from those In 
the past studies In that convection Is considered to be of primary Importance 
In the present model. Consequently, the coupling effect mentioned above can 
no longer be Ignored. Among several Important and Interesting features of th* 
thermocaplllary flows, the surface temperature distribution Is one of the most 
critical factors because It determines the driving force of the thermocaplllary 
flows. In contrast to the conduction-dominant situations, the surface temper- 
ature distribution due to thermocaplllary convection In a thin liquid layer 
with a nonuniform heat flux Imposed on the free surface Is determined by the 
coupling among the flow field, the temperature field and the Imposed heat flux. 
The surface temperature distribution Itself Is a part of the solution. Con- 
sequently, any estimation of the global surface-temperature gradient based on 
the Information derived from the conduction-dominant situations will be Inap- 
propriate. Without the correct Information about the global surface tempera- 
ture distribution, a proper estimation of the driving force and, hence, other 
Information about the thermocaplllary flow cannot be obtained. Furthermore, 
for a two-dimensional thermocaplllary convection experiment. If the dimension 
In the flow direction Is not designed large enough, sharp surface temperature 
variations In the corner regions near the lower temperature side walls will 
exist. This sharp surface temperature variation will bring about strong shear 
stresses In corner regions, which, under certain circumstances, may Induce 
unsteady flows. Lowry (1980) experimentally studied a similar problem using a 
two-dimensional container with a horizontal heating element located symmetric- 
ally above the free surface of the fluid sample. The velocity and temperature 
distributions were measured for several different cases, but the relationship 
between the surface temperature distribution and the coupling was not taken 
Into consideration. As a result there has been no guidance for modeling or 
designing an experiment to avoid sharp surface-temperature variations which are 
known to exist In the corner regions near the lower temperature side walls. 

The present study Is aimed at resolving all these deficiencies. 

Dimensional analysis with proper balancing between physical quantities Is 
employed Initially to obtain appropriate dimensionless parameters. Based on 
these parameters, the characteristic temperature difference along the free 
surface AT 0 and the characteristic length scale In the flow direction L 
can be determined. Having had these two quantities defined, the global 
surface-temperature gradient AT 0 /L, the global driving force, and other 
Interesting features such as the characteristic velocity can be realistically 
assessed. A numerical calculation Is used to verify the validity of the param- 
eters obtained In the first part of the present study. Hopefully, this effort 
will lead to a better understanding of certain Important and Interesting fea- 
tures of thermocaplllary flows. 
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MATHEMATICAL FORMULATION AND DIMENSIONAL ANALYSIS 


As a result of the coupling among the flow field, the temperature field, 
and the Imposed heat flux, the surface temperature becomes a part of the solu- 
tion Itself when convection becomes Important. The purpose of this section Is 
to obtain a proper estimation of such a surface temperature distribution. 

The schematic diagram of the configuration and the applied coordinate 
system of the present study are shown In figure 1. The problem considered 
herein Is a steady, two-dimensional thermocaplllary convection In a thin liquid 
layer (l.e., Inertia effects are negligible) with nonuniform heat fluxes 
Imposed on the free surface. For simplicity, but without losing the general- 
ity, a Gausslan-dlstrlbutlve heat flux Is considered In the present study. The 
bottom boundary Is kept at the lower ambient temperature. Under the thin layer 
configuration, buoyancy force Is negligibly small compared to the thermocapll- 
lary driving-force. The dimensionless governing equations and boundary condi- 
tions are as follows: 
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where (u,v) Is the dimensionless velocity vector, p Is the dimensionless 
pressure, e Is the dimensionless temperature difference, A Is the aspect 
ratio, and Ma Is the Marangonl number. In the above formulation, the free 
surface Is assumed to be flat, and the heat loss from the free surface Is 
neglected. In addition, the following nondlmenslonallzatlon scheme has been 
employed: 


x ~ L , y ~ D, u ~ — 
k l r 


0 |j t 1at r 

V 


o£ |q T |AT R 

~ L 2 v 

L R 


P ~ 


| °T | AT R 


0 ~ AT, 


3 



where y denotes the dynamic viscosity, oj Is the surface-tension gradient 
with respect to temperature T, D Is the layer depth, Lr Is the character- 
istic length scale of the Imposed heat flux. Let f 0 be the amplitude of 
the Gausslan-dlstrlbutlve heat flux and k be the thermal conductivity. Then 
ATr - f o 0/k gives the referenced temperature difference. It Is derived 
from the energy boundary condition at the free-surface by assuming conduction 
as dominant. 

Ti.e dimensionless parameters Involved are 

aspect ratio A « p- 
l R 

OkjhT 

Harangonl number Ma * - 

ya 

where a Is the thermal dlffuslvlty. When convection becomes Important, that 
Is, when (A 2 Ha) » 1, Lr and ATr are no longer the proper characteristic 
quantities because they are derived under conduction dominant conditions. The 
global surface-temperature gradient based on the ratio of these two quantities, 
ATr/Lr, becomes therefore Incorrect. To obtain appropriate estimations of 
a proper length scale L and a proper characteristic temperature difference 
AT 0 In terms of Lr, ATr, and A 2 Ma, the following argument Is made. 

In the steady state, L and AT 0 can be properly determined from the 
global energy-balance condition (l.e., the energy Input from the Imposed heat 
flux should be equated to the energy output by conduction through the Dottom 
wall) and the local energy boundary condition at the center point (0,1) of the 
free surface, where che temperature reaches Its maximum value. The global 
energy-balance condition Indicates that 


or 





(f 0 D/k) 



(9) 


The expression k(AT 0 /D) has been used above to estimate the global heat flux 
through the bottom wall. This Implies that. In the steady state, convection 
should finally be balanced by conduction across the entire layer. Otherwise, 
convection will keep heating up the liquid In the downstream until the flow 
field eventually reaches the steady state. But at the central position, where 
the thermal boundary layer has not been developed to the full layer depth, the 
use of D as the characteristic length scale for the conduction In 
y-dlrectlon will be Improper. The local energy boundary condition at the 
center point of the free surface (eq. (6)) Indicates that 
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where denotes the thermal boundary layer thickness at the central pos- 
ition, which Is obtained by balancing conduction and convection In the 
y-dlrectlon at that position. Let 
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be the characteristic velocities In the x- and y- directions, respectively, 
when convection becomes Important, and A* - D/L and Ha* > (lljl/a), the 
corresponding modified aspect ratio and Harangonl number; then, by the balance 
between convection and conduction In the y-dlrectlon at the central position, 
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and (using eq. (9)) 
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and equation (10) applies. 
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Now we have the following two relations: 
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and 


»T 0 - <A 2 H»r V4 »T R (12) 

Equatlo > (11) and (12) Indicate that L Is directly proportional to (A 2 Ha)^* 
and aT 0 Is Inversely proportional to (A 2 Ha) 1/4 when convection Is Important. 
These two relations also Indicate that as convection becomes more and more 
dominant, that Is, when the value of A 2 Ma becomes larger, more heat gets 
transported downstream by thermocaplllary convection. As a result of this, the 
characteristic temperature difference Is reduced to a smaller portion of ATr. 
These two relations are shown In figure 2 and will be verified In the next 
section when numerical calculations are presented. 

A further remark should be made with regard to the values of A* 2 Ma*. 
Combining equation (10) with equation (11) or (12) yields 

A* 2 Ma* ~ (A 2 Ma) 1/4 


This Implies that, In the steady-state, although the convection should finally 
be balanced by the conduction In a thin liquid layer as discussed earlier. It 
Is Inadequate to require A* 2 Ma* ~ 1. The unusual result Is due to the 
deceleration of the surface velocity In thermocaplllary flows ( 1 . e . , from 
Uq to zero at a far downstream point of the flow field), which allows the 
thermal boundary layer to grow faster than expected. 


NUMERICAL CALCULATIONS 


In the thin liquid layer analysis, as the velocity field Is coupled with 
the temperature field only through the shear-stress boundary condition at the 
free surface (eq. (5)), the momentum equation (eq. (2)) can be solved first. 
The solutions of u(x,y) and v(x,y) are 


u(x,y) 



3 2 



and 


v(x,y) 





where e s (x) Is the dimensionless surface-temperature distribution. The 
expressions of u and v also satisfy the continuity equation. The energy 
equation (eq (4)) becomes, after substituting the above expressions for u 
and v. 
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Because of the local nonlinearity (l.o., the existence of [(de s /dx)(a«/»x)] 
and [(d 2 e s /dx 2 )(d6/dy)] In the convection terms). It Is still very difficult 
to obtain the analytical solution of e(x,y). In this section, a finite dif- 
ference numerical scheme Is employed to calculate the temperature distributions 
•(x,y) for different values of A 2 Ma. 

Instead of solving the steady state equation, a time dependent term Is 
Included for the time marching procedure. The equation to be solved becomes 


»? * * 2n * 


de 




dx 


d 2 e 


z * 1 « * 4 (y3 - y2) If 


& 

ay 2 


The main features of the numerical scheme employed In the present study can be 
outlined as follows: 


(1) Central differences are used for both the first and second derivatives 
except at the boundaries where suitable boundary conditions are Incorporated. 

(2) Because of the modest variations, the grid spacing for the numerical 
calculation are chosen to be 0.1 In both xand y-dlrectlons. 

(3) Depending on the magnitude of A 2 Ma, proper time Intervals are 
chosen for the time marching procedures so that stable numerical calculations 
can be obtained. 


(4) To determine the steady-state, a proper criterion c Is used for the 
temperature difference calculated at two time steps. The criterion depends on 
the magnitude of A 2 Na and the global energy balance between the Imposed 
heat flux and the heat conduction through the bottom wall (to be within ~10 -2 
difference). As an example, the criterion c used for the case with A 2 Ha 

■ 10* Is ~10~ 3 for two calculations at one-thousand time stejs apart. 

(5) In order to reduce the computation time, the domain of calculation Is 
Increased successively when appropriate. That Is, when the temperature at the 
last grid point reaches 0.002, another ten more grid points at the downstream 
are then added to the calculation domain until the temperature field finally 
reaches the steady-state. The Initial domain used for calculation Is 30 grid 
spaces In the x-dlrectlon (l.e., up to x - 3) and 10 grid spaces In the 
y-dlrectlon (l.e., up to y ■ 1, the free surface). 

Since the detailed temperature distributions of the whole flow field are 
not of Interest to the present study, only the surface temperature distribu- 
tions for various situations are recorded and plotted In order to verify the 
theoretical predictions made In the previous section. The surface temperature 
distributions for various values of A 2 Ma are shown In figure 3. It can be 
seen that as the convection become more Important, that Is, as the value of 
A 2 Na Increases, the extent of the surface temperature distribution becomes 
larger. Consequently, the dimensionless temperature difference at the central 
position decreases, that Is, the characteristic temperature difference AT 0 
becomes a smaller portion of the referenced temperaii're difference aTr. 

These two behaviors had been predicted by nond Intensions! analysis discussed In 
the previous section. The relations of the characteristic length scale L 
and the characteristic temperature difference AT 0 versus A 2 Ma are 
shown In logarithmic scales In figure 2. Although the numerical calculations 
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of the characteristic temperature difference do not appear to coincide with the 
theoretical line, they are within the same order of magnitude. In the non- 
dimensional analysis with order-of-magnltude estimates, the power law Is of 
primary concern. It can be seen from figure 2 that both the theoretical 
predictions and the numerical calculations follow the same trends. This 
confirms the validity of equations (11) and (12). 


SUMMARY AND CONCLUSIONS 

The thermocaplllary convection in a liquid phase with a nonuniform heat 
flux Imposed on the free surface Is a vt r v practical and Important problem 
related to the floating-zone c;ystal growth. For simplicity, a two-dimensional 
thin liquid layer Is employed herein to Investigate the effect on surface 
temperature distribution due to the coupling among the flow field, the temper- 
ature field and the Imposed heat flux. The following conclusions are drawn 
from the analysis made In the present study: 

1. Dimensional analysis and the numerical calculations Indicate that for 
thermocaplllary convection flow In a thin layer configuration, the following 
two relations are valid: 


L 


(aW / 4 l, 


and 


2 -1/4 

AT q ~ ( A*Ma) ' AT r 

With L and AT 0 determined, the global surface-temperature gradient, the 
thermocaplllary driving force and other Important Information such as the 
characteristic velocity can be assessed properly. 

2. With the above two relations, proper experimental designs for the 
Investigation of the thermocaplllary convections can be made to avoid the sharp 
surface-temperature variations In the corner regions near side walls. Effects 
due to strong shear stresses In corner regions near the side wills can then be 
reasonably eliminated In such a configuration and If an unsteady, oscillatory 
motion Is found. It can be attributed to coupling among t'.e flow field, the 
temperature field and the Imposed neat flux. 
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